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Abstract 

Ab initio atomic relaxations often take large numbers of steps and long times to converge. An 
atomic relaxation method based on on-the-flight force learning and a corresponding new curved 
line minimization algorithm is presented to dramatically accelerate this process. Results for metal 
clusters demonstrate the significant speedup of this method compared with conventional conjugate- 
gradient method. 
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One major usage of ab initio density functional theory (DFT) in material science simula¬ 


tion is to determine the ground state atomic conhguration for a given system [i|, [^. Overall, 
such applications probably take most of the DFT simulation time. There are two types of 
ground structure searching. The hrst is to hnd global minimum among many local minima 


j3|, 1^. This has become an intensely studied topic in material desi gn p rojects [5l-l9|. Vari¬ 
ous types of evolutionary algorithms jb-lO] or simulated annealing 11 1 schemes have been 


developed, as well as the minimum hopping methods The second type is the conven¬ 
tional local minimum optimization, which is the concern of the current study. The related 
calculation is dominated by the conjugated gradient (CG) method jl2l-ll4| and the Broyden- 
Fletcher-Goldfarb-Shanno (BFGS) method 15l-ll7l|. Although these methods guarantee to 
converge into a local minimum, the convergence rate could be agonizingly slow, e.g. with 
hundreds of steps, thus a faster method will be extremely helpful. This local minimum prob¬ 
lem also presents itself in the global minimum search since each global minimum search step 
usually deploys one or more local minimizations jd-O]. One reason for the slow convergence 
of the local minimization steps is the possible narrow and curved energy valley leading to 
the minimum, which prevents the efficient execution of the conventional GG or BFGS meth¬ 
ods. Imaging a rotation of a molecule on the surface of a substrate. Such rotation cannot 
be described by a straight line in cartesian coordinates which is used under GG or BFGS 
methods. In higher dimension, the situation can be more complicated, making it impossible 
to hnd the natural degree of freedom (e.g., the rotation angle). One such example is a metal 


cluster 


18 


19| (which will be studied in this paper), where hundreds of steps might be 


needed to relax a structure while there is no obvious natural (or say internal) degree of 
freedom to speed up the convergence. To overcome these problems, one needs to do the 
minimization steps along guided curved lines following the energy valleys. We will call such 
algorithms the guided curved-line-search (GLS) algorithms. 

The issue is how to hnd such guided curved lines. In this work, we will show that such 
guided curved line can be provided by model surrogate potentials with their parameters 
provided by on-the-hight htting (OTFF) to the ab initio atomic forces 20l-l22| . We will 
demonstrate the efficiency of our GLS algorithm on metal clusters. Overall, we have the 
following hndings: (1) The GLS method can speed up the traditional GG method by a factor 
of 3 to 6 for both the number of steps and wall clock times; (2) The OTFF can be ehectively 
used to speed up the atomic relaxations, not just molecular dynamics as it has been used 
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so far 
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22I; 


As mentioned before, the guided curved line will be provided by a surrogate potential. 
One possible option is to carry out ab initio line minimization along the steepest descent line 
(SDL) of this surrogate potential. We will use on-the-flight htting (OTFF) to ensure that 
the atomic forces of this surrogate potential at the beginning of each step equal that of the 
ab initio forces. When the system approaching the hnal minimum point, the curved line will 
become straight in the small scale, then the curved line search will go back to the conventional 
straight line search. In practice, we found that the SDL can be warped with sharp twists 
in high dimensions. Besides, using SDL will miss the conjugated gradient feature between 
different line searches. To overcome these shortcomings, we will use the surrogate potential 
conjugate gradient descent line (SP-CGDL). To construct SP-CGDL, the conventional GG 
formalism is applied to the initial atomic force direction to yield the GG search direction. 
Then a straight line minimum search based on the surrogate potential is carried out. From 
the new line minimum point of the current surrogate function, subsequent GG straight lines 
are carried out. Thus, our SP-GGDL curved search line is consisted with many straight 
lines segments. The ab initio line minimization will be carried out alone this SP-GGDL. 
One might worry that the ab initio energy function along this segmented line might not be 
smooth enough to carry out ab initio line minimization. But in practice, we found that one 
can effectively use the Brents algorithm 231] to search for the ab initio line minimum along 
this SP-GGDL, and such line search often hnds the line minimum at a few segments down 


the road along the SP-GGDL. Typically two ab initio calcu 
algorithm to search for the minimum along the SP-GGDL 


ations are needed in the Brents 


24l |. much like the conventional 


line minimization calculation. After the ab initio line minimization are done, we call this 
one step, and the algorithm will repeat itself (from OTFF to construction of SP-GGDL, 
then ab initio line minimization). Note, the above procedure maintains the feature of GG, 
when close to the minimum. 

We will use metal cluster 
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19| to demonstrate our GLS algorithm. The metal cluster 


potential is intrinsically high dimensional due to their long range atom-atom interaction. As 
a result, it is often difficult to reach their local minima. The metal cluster is an important 

18| . A lot of works have been done in searchin g of the 
optimal cluster structures, and the density of local minima in energy 


subheld related to catalysts 
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For metal systems, we found that the Wbody Gupta force held 


28| is a very good general 
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potential. It has been used to model various types of metal clusters. The potential is a special 


case of the embedded atom potential 


29l | based on the second moment approximation of the 


tight binding theory and it has the following form: 
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where represents the distance between the atom i and j in the cluster. The hve groups of 
parameters Aij, ^ij, pij, qij, r°- are allowed to vary independently to match the forces from 
DFT calculations. We restrict Aij = Aji, ^ij = pij = pji, q^j = qji and r^j = r^i, and 
set them to zero when the is larger than a cut-off distance to limit the number of the 
variables. Parameters are also restricted to vary within a physically meaningful range. 

Because the analytic expression for atomic forces of this model is a non-linear function 
of these parameters, in order to have an accurate force htting, we have used a parallel 
differential evolutional algorithm 10| to globally minimize the force error. The resulting best 
solution is further optimized by a CG local minimization algorithm for these parameters. 
This approach enable us to always £t the atomic forces with an error less than 0.005 eV/A, 
which is a few times lower than the typical ab initio minimization stoping criterion. Although 
the fitting procedure (at the beginning of every ab initio line minimization step) might sound 
complicated, its computational cost is negligible, about 5% of the ab initio computational 
time. 

To show the quality of atomic force htting, we present the atomic force error in FigH] 


for a P 
Ref. 281 




cluster (with 100 Pt atoms). To begin with, we use the Gupta parameters from 


30l . and l3l| which have the parameters for almost all the major metallic elements. 


The atomic force error compared to ab initio calculation using these original parameters 
without htting is about 1 eV/A. After the parameter htting, they becomes about 10“^ 
eV/A. This improvement on the force is at no cost of degradation of other properties of this 
functional. For example, FiglT](c),(d) compare the atomic force changes between the Gupta 
and DFT results when the atomic positions have been randomly displaced. The original 
Gupta result is already rather good, and it has been slightly improved after force fitting^ 


To demonstrate the speedup of the GLS method, we hrst test hve random Pt 2 o clusters [10 1 
with diherent initial structures and corresponding diherent initial energy. The convergence 
results are shown in Figj2](a) in comparison with the conventional GG results. We see that. 
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FIG. 1. Deviation of atomic forces for Ptioo between DFT calculations and approximate models: 
(a) conventional un-fitted Gupta potential; (b) force-fitted Gupta potential. Comparisons of the 
values of force difference F — Fq between DFT and approximate models: (c) un-fitted potential; (d) 
force-fitted potential, where Fq is the force at the atomic structures (i?o) used for force fitting and 
F is the one at randomly displaced structures {R) around Rq. Color bar for (c) and (d) indicates 
the distance (in unit of A) of R from Rq. Note that the plot contains several different groups of 
Rq. In (a)-(b), the x-axis is DFT force while the y-axis is the deviation (in unit of eV/A). In (c) 
and (d), the x-axis is DFT force difference F — Fq while the y-axis is the approximate one (in unit 
of eV/A). 

more than a factor of 3 speedup are achieved for most cases. Especially, more than a factor 
of 6 speedup are achieved in the initial relaxation steps. Note, the total time for each line 
minimization step for CLS and the CG method is approximately the same. 

In actual work, one often uses Gupta to pre-relax the system to a Gupta local minimum, 
then uses conventional ab initio GG relaxation to further relax the total energy of the system. 
We will call such scheme pre-GG. One can also start with the Gupta relaxed minimum, 
then use our GLS method, we will call such method pre-GLS. Their results are shown in 
Figj2](b),(c). We can see that, due to the good approximation of Gupta to DFT energy, 
the Gupta minimum does provide a good initial speedup for ab initio energy minimization, 
although pre-GLS still out performs pre-GG by a factor of 2 to 4, and if accurate results 
are needed, subsequent minimizations are as costly as the original GLS and GG methods. 
For the global minimization search problems, or to search for local minima density, one 
issue is that the pre-GG or pre-GLS method tends to mislead the system to the same local 
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FIG. 2. (a) The relaxation process of five random Pt 2 o clusters with different initial structures and 
corresponding different initial energies by the CLS and CG method, respectively; (b) the convergent 
process of the mean energy of these Pt 2 o geometries after pre-relaxation, and the relaxation of Ptioo, 
(c) Coi 2 o and Cu 2 oAui 8 clusters with the initial geometry from the global minimum of conventional 
un-fitted Gupta force field. The x-axis is the number of relaxation steps while the y-axis is E — Ef 
(in unit of eV), where E is the energy of the current step and Ef is the energy of the finally sought 
structure. 


minimum near the Gupta potential basins for different initial configurations. This reduces 
the diversity of the global search. On the other hand, we found that the CLS method is 
exempted from such a problem. 

Finally, in Figj2](c), we show that the CLS method also works for Co and CuAu alloy 
clusters, demonstrating its generality for metallic systems. 

In summary, we have presented a curved line search (CLS) algorithm to speed up ab initio 
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atomic structure relaxation. This CLS uses a classical potential to provide the curved line 
on which ab initio line minimization is carried out. The parameters of this classical potential 
are htted on-the-flight at every step to the ab initio atomic forces. We tested this approach 
using metal clusters with Gupta force held as the classical potential and we expect similar 
approaches can be applied to other systems. Compared to the traditional CG method, we 
found CLS can speed up by a factor of 3-6. 
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